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Abstract. H2O is a key ingredient in many interstellar environments, like photon dominated regions and star forming clouds. It 
plays an important role in the oxygen chemistry and can act as a coolant in dense interstellar clouds and shocks. Observations 
and modelling of water lines thus provide powerful diagnostics of the physical conditions in interstellar emission zones. A 
radiative transfer method for the treatment of molecular lines is presented. We apply this method to previous SWAS and ISO 
observations of water vapor in the source S 140 in order to make models to plan for, and to interpret, HIFI data. Level populations 
are calculated with the use of a three-dimensional (multi-zone) escape probability method and with a long characteristic code 
that uses Monte Carlo techniques with fixed directions. Homogeneous and inhomogeneous models are used to compute the 
differences between HtO line profiles across the S 140 region. We find that when an outflow or infall velocity field with a gradient 
of a few km s"' is adopted, line profiles with a FWHM of 6 km s"' are found, in agreement with observations. Inhomogeneous 
models are favoured to produce a single-peaked line profile. When zooming in on smaller regions within the PDR, the shapes of 
the line profiles start to differ due to the different temperature and density distributions there. The embedded sources are traced 
by high excitation lines of, e.g., 32i - 221, 3o3 - 2i2, 2,2 - loi and 220 - lii. The computed intensities are roughly consistent with 
existing ISO observations. Water emission in a PDR source like S140 requires a combination of a pure PDR and an embedded 
source in order to match the observations. Because of its good angular resolution, HIFI will be able to distinguish between a 
dense star forming region or a more diffuse gas component. It is therefore important for future observing programs to consider 
both in their predictions of the emission characteristics of water in these environments. 
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•I— I 1. Introduction wing) could dissapear when mol ecules are deplete d in spite of 

infall velocity fields. Moreover, ' Lee et all ll2004h argued that 

Our understanding of astrophysical objects, such as molec- ,• , 1. ■ j . j r j ^ j 

3 , , , , . .• . , line observations can be misunderstood if one does not consider 

^« ular clouds, young stellar obiects (YSOs), disks, plane- , • • ■ 1 j ■ , ■ • 1 j- xt • 

, , , , , • , • /T,.^T, ^ • the variations in abundance with time or with radius. New in- 

tary nebulae and photon dominated regions (PDRs), requires ^ ^ j j , • • j n 

J- J . . J. . struments are needed to zoom in on these regions and spatially 



both observations and, in order to retrieve information on 



resolve them in order to retrieve more accurate information on 



physical parameters, detailed modelling. Line emission from • 1. • 1 j 1. 1 j-.- t-u 1 u cn u 1 

^ ^ their physical and chemical conditions. The launch of Herschel, 

molecules, such as CO and H2O, has been observed exten- 1 j rr . j t r i t r j 

J J l„ . — — I 11— — ■ — 7— -1 with on board the Heterodyne Instrument for the tar Infrared 



I ' 1 -ii — 1 witn on Doara tne neteroayne instrument for tne tar infrarea 

sively in the past decades JFederman et all 19801 Helmich et alJ /titt-t\ n 1 ui . u .u i- t f j u 

. '~J 17- : „ ^ . , ' L^^J , .. ,M,„„„I (HIFI), will make it possible to observe the Far Infrared world 

in more detail than ever before with higher angular resolution 



1996; Herpin & Cernicharo 2000; Hogerheiide et al. 1998' 
Ladd et al. 1998; Stark et al. 2004; van Dishoeck 2004) and has 
been proven to be a powerful tool to reveal the behaviour of the 



and sensitivity. As more infrared and submillimeter surveys are 

planned in the future (e.g., ALMA), more models to interpret 
aformentioned environments. Line profiles of molecules con- ^, , ^ c ■ ^ jj 

. .. . „ . the data from these instruments are needed. 



tain important information on the conditions of the regions 
they occupy, such as the density, and the velocity structure 



One way to model the line emission from various molecules 

^ . . . ^ . ^ , ,. . , is to make use of an escape probability method, in which the 

Retrieving information from these lines requires accurate mod- t. t.-t^ r, ■ . i, c r. . ^ .i. 

„. , . . r , , . , probabihty a indicates the chance for a photon to escape the 

elling. However, the interpretation of the line shapes has to be i j rr-,- j • u- u i i n .i. 

^ 1 J r ■ cloud. This method, in which one solves locally the equa- 

treated with great caution. lRawlings & YatesU2001D noted that . ^ . , ■ j • . u u j 

" ^ ' ^ tions of statistical equilibrium in every gridpoint, has been used 

the classic signature of infall (an asymmetric double-peaked • . . i i . .v. j- i- r 

, . , , , . , . , , , in the past to calculate, e.g., the radiati ve coohng of warm 

line which possesses a blue wing that is sti'onger than the red , , jrr — n „ c linmlv .u r • r j 

° ° molecular gas (Neu feld & Kau fman 1993i. the far-infrared wa- 

Send offprint requests to: D.R.Poelman® astro.rug.nl ter emission from shock waves (.Kaufman & Neufeld 1996.) and 
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water maser emission jNeufeld & Melnicl3ll99ll) . Despite the 
wide range of applications of this method, other techniques, 
e.g., Monte Carlo (MC), Accelerated Lambda Iteration (ALI) 
and Local linearization (MULTI type) codes, that solve non- 
LTE multi-level radiative transfer problems, were created in 
order to meet with the non-local effects that arise in certain 
problems. Monte Carlo codes have been used over the years 
to study, e.g., masers ( Spaans & van Lanaevelde 1992), UV 
continuum transfe r (Scaans 1996) and molecular line transfer 
llHogerheiide & van der Tak,2Q0Q) . These have been proven to 
be reliable tools. However they have the disadvantage of hav- 
ing conve rgence problems at high op tical depths (t > 100). 
Recently, lElitzur & Asensio RamosI ( I2006i) introduced a new 
method to treat the radiative transfer of lines. This method, 
the Coupled Escape Probability (CEP), combines the algebraic 
equations used in the escape probability approximation with a 
closed set of equations in which the effect on the level popu- 
lations in some zone / of radiation produced in all other zones 
is taken into account. Hence, the CEP method retains all the 
advantages of the escape p robabihty approach. An e l oquen t 
discussion can be found in Elitzur & Asensio RamosI (|2006). 
Our method, developed independently of the CEP, resembles 
the latter in that we also make use of a multi-zone formalism. 
However, the non-local component i n the CEP treatment, i.e., 
the sec ond term of equation (35) in Elitzur & Asensio Ramoa 
is replaced in our formalism by a spatial iteration loop 
on the different escape probabilit ies in ever y zone along a num- 
ber of discrete directions ( Poelman & Spaans..200S,) . Note as 
well that our formalism already is fully 3D. 
In this work, the multi-z one escape probability method of 
IPoelman & SpaansI ll2005l) is used to compute, through a ray- 
tr ace technique, line profiles for the S140 region. The models 
of |Poelman & SpaansI I2OO?) for the S140 PDR are extended 
through the inclusion of embedded sources. For details on the 
physical structure of the S140 P DR and the employed escap e 
probability formalism, we refer to P oelman & S paans (2005). 
In Section |2| we discuss the excitation of certain water lines 
in an embedded source. Section |3 describes the line profiles 
in S140 for the ground state transition of ortho-water. Also, in 
the appendix we test our escape probability model against a 
benchmark Une profile problem and present a long characteris- 
tic radiative transfer method. 

2. Embedded Sources 




ID'S 10'' 10'= 10" 

Distance from central star (cm) Distance from central sfcn (cm) 

Fig. 1. The distribution of density (left) and temperature (right) 
as a function of position from the embedded source. 




10"s 10" lO's 10" 

Distance from central star (cm) Distance from central star (cm) 



Fig. 2. The water abundance distribution, relative to H2, as a 
function of position from the embedded source for two (model I 
(left); model II (right)) different scenarios. 

Table 1. Observed ISO line intensities [erg cm"^ s ' sr '] for 
S140 IRSl 



Transition and Wavelength 

113.5 yum 136.5 yum 174.6 yum 179.5 yum 

4[4 — » 3o3 330 — > 321 3o3 — » 2i2 2[2 — » Ipi 

<6.8xlO-« < 8.6x10-" 7.1x10-'' < 4.5 x lO-"* 



Before going to Sect. |3 in which we discuss the line pro- 
files for the S140 region, we calculate the radiative transfer 
of certain transitions of ortho- and para-H20 when an embed- 
ded source i s incorporated into t he mo dels of the S140 PDR 
presented in IPoelman & SpaansI 120051) . The PDR lies at the 
south-western edge of the molecular cloud LI 204. The molec- 
ular cloud, which extends over more than 30' (10 pc), lies at a 
distance ~ 900 pc, and is illuminated by a BOV star HD2 1 1 880, 
located at about ~ T to the cloud. As known from observa- 
tions (IRS 1 to 3. lEvans et aL. 1989) . an embedded cluster of 
infrared sources is part of the S140 molecular cloud. The tem- 
perature and density structure characterizin g these sources are 
different from the ones adopted in Poelman & Spaan^ ( l2005l) 
for the PDR and the Extended Molecular Cloud (EMC), and 
is therefore a useful extension to the performed calculations in 
the aforementioned paper. 

The calculations have been performed in ID (201 gridpoints), 
with an adopted size of 0.06 pc and 0.12 pc for the em- 
bedded source in model I and II, respectively. Using spher- 
ical symmetry, the resulting level populations have been ro- 
tated to become a 3D solution from which the intensities have 
been computed. No systematic velocity field is included, but 
a velocity dispersion, Aves = (\h ^urb-*'^^' ^ km s"' is 
adopted. The temperature and density structures of the embed- 
ded source were derived bv . van der Tak et al. ( 2000. 2003) , us- 
ing a power-law density structure n = «o(r/ro)~". The adopted 
density and temperature profiles are plotted in Fig.[l] Following 
[ Soonman et al. ( 200 3), our density increases from a few times 
10^ cm~^ at the edge to almost 10^ cm"^ in the centre of the 
embedded source for both models. We assume, for simplicity, 
that the temperature of the gas is coupled to the dust because 
of the high densities. In model I and II, the temperature rises 
from, respectively, 80 K and 60 K at the edge towards 450 K 
in the centre. The corresponding H2O abundance profiles are 
plotted in Fig.|2| Note that the average H2O abundance in the 
PDR is ~ 1 X 10"**, whereas it is a factor 20 and 2 x 10^ higher 
in Model I and II, respectively. The adopted water abundance 
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transition 




Wavelength 


average 


average 


average 








intensity 


intensity 


intensity 








model I 


model II 


Poelman & Spaans (2005), FDR I 






\pLm\ 


r -I -In 

ere cm " s sr 


Ferg cm"^ s"' sr"'l 


r -1 -In 

em cm " s sr 




loi 


538.3 


5.0 X 10"^ 


2.2 X 10"* 


1.1 X 10"' 


O-H2O 2i2 -> 


loi 


179.5 


4.5 X 10"' 


1.4 X 10"' 


6.2 X 10"^ 


0-H2O 221 -» 


1,0 


108.0 


3.4 X 10"' 


1.4 X 10"'' 


5.7 X 10"' 


0-H2O 221 -» 


2,2 


180.5 


1.1 X 10"' 


3.8 X 10"' 


4.1 X 10"'" 


O-H2O 3o3 ^ 


2l2 


174.6 


2.5 X 10"' 


5.4 X I0-' 


4.7 X 10"' 


O-H2O 3i2 -» 


221 


259.9 


6.8 X 10"* 


3.1 X 10"' 


1.3 X 10"" 


O-H2O 3l2 -» 




273.2 


4.9 X lO"* 


1.6 X 10"' 


7.8 X 10"" 


O-H2O 321 


2l2 


75.3 


3.7 X 10"' 


2.2 X 10"'' 




0-H9O 391 -> 


3l2 


257.7 


2.7 X 10"* 


8.0 X 10"* 




0-H2O 4i4 -> 




113.5 


2.7 X 10"' 


1.3 X 10"" 




0-H2O 3to -» 


321 


136.5 


5.0 X 10"* 


2.2 X 10"* 




P-H2O ill -» 


Ooo 


269.3 


1.5 X 10"' 


1.5 X 10"' 


2.5 X 10"' 


P-H2O 2o2 -» 


In 


303.4 


1.2 X 10"' 


1.0 X 10"' 


2.7 X 10"** 


P-H2O 220 


In 


100.9 


1.7 X 10"' 


1.4 X 10"" 


8.8 X 10"'" 


P-H2O 2i, ^ 


2o2 


398.6 


4.2 X 10"* 


4.6 X 10 * 


2.7 X 10"'" 


P-H2O 220 ^ 


2ii 


243.9 


2.2 X 10"* 


1.4 X 10"' 




P-H2O 3 13 


2o2 


138.5 


1.8 X 10"' 


7.8 X I0-' 




P-H2O 3 13 ^ 


^20 


1635.7 


1.9 X 10"" 


1.6 X 10"^ 




P-H2O 322 ^ 


2,1 


89.9 


1.8 X 10"' 


1.5 X 10"" 




P-H2O 4o4 


3,3 


125.3 


1.4 X 10"' 


9.5 X 10"' 





^ Like^o^^^Jeufeld ^9T) we find weak maser action in this line. 



profiles are merely to illustrate the effects that are expected to 
result when wat er molecules are thermally released from warm 
dust grains (c.f. iBoonman et aDcOOS) . Note that our embed- 
ded source s are only ~ 2 X lO'^ cm in size, compared to ~ 
10'^ cm in lBoonman et alJ 12003). We feel that this is appro- 
priate becai ise 10^^ cm is approx imately the size of a clump 
in the PDR jPoelman & Spaansl2005b . He nce, beyond lO'^ cm 
photo- dissociation of water as described in Poelman & Spaanj 
( I2OO5I) dominates. Table [T] lists th e obser ved ISO line inten- 
sities discussed in iBoonman et alJ ll2003h . The observed line 
fluxes [W cm"^ Afm"^] have been converted into line intensi- 
ties [erg cm"^ s"' sr"^] by adopting a beam size of ~ 78" (~ 
0.4 pc at the distance of S140). The resulting intensities of the 
two adopted models are listed in Table |2] together with the re- 
sults found in Poelman & Spaans (2005) for the PDR without 
embedded sources. Higher densities and temperatures lead to 
increased intensities for high lying transitions. The lio — > loi 
at 557 GHz is a factor ~ 10 too high compared to the observed 
SWAS value for both the mod els. However, beam dilution ef - 
fects are likely to alleviate this (' Spaans & van D ishoeck'200l'). 
The predicted intensity for the 3o3 2i2 line at 174.6 jum is 
less than a factor of two lower than the observed value of 7.1 
X 10"^ erg cm"^ s"' sr"' for model II, while model I is lower 
by a factor of ~ 3. A core of water, as in model II, is preferred 
to model this transition. The predicted intensities for the 2,2 — > 
loi line at 179.5 fim and the 330 ^21 line at 136.5 yum are in 
good agreement with the observed upper limits for both mod- 
els. The observed upper limit for the 4i4 — > 3o3 line at 113.5 
//m is best reproduced by model I. Model II produces too much 
intensity for this transition by a factor of ~ 20. 
The intensities listed in Table|2]show that water emission from 



a PDR source, like S 140, requires a combination of a pure PDR 
contribution for 0-H2O lio — » loi dPoelman & S paans"2005|). 
and an embedded source for high excitation lines (Model I/II) 
in order to match the observations. 

3. Line profiles 




Fig. 3. The figures show the predicted distribution of the intensi- 
ties (erg cm"^ s"' sr"') of two transitions of ortho-H20. The left 
figure depicts lio — > loi, on the right one has the 3,2 — > 221 
transition. The star is located to the left of the figures. 

In this section, we present line profiles for the lio — > loi 
transition of ortho-H20 at a few positions in the S140 PDR 
(without embedded sources), as well as line profiles when the 
beam is as large as the PDR (Fig.|5]&0. Line profiles for a 
static embedded source (Fig.|6j are also shown. 
The computed temperature, density and water abundance dis- 
tributions, as calculated with the inhomogeneous code of 
Spaans (1996), is used as input for solving the H2O equa- 
tions of statistical equilibrium by means of an escape proba- 
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Table 3. PDR model parameters 



Model 


size 




n- ^ 

"■iC 






Av/ 




[pc] 


[cm-3] 


[cm-3] 


[%] 


[pc] 


[km s-'] 


PDR 


0.5 


4x 10^ 


1 X 10* 


2.5 


0.03 


1.2 



" Clump density 
* Interclump density 

Fraction of the total volume that is occupied by the clumps 

Clump size 
" Velocity dispersion 




Velocity (cm s->) 



Fig. 6. Line profiles for the lio — > loi (dashed) and 3o3 — > 
2i2 (solid) transitions of 0-H2O for a static embedded source 
(Model 1). Y-axis in units of erg cm"^ s"' sr"' per cm s"'. 

bility method. For a full d escription of the model, we refer to 
IPoelman & Sp aans (2005). Fig.|3lshows illustrative maps of the 
predicted distribution of the intensities of two ortho-HaO tran- 
sition s^_TTie_2^atnetersjised^ this model are listed in Table 
l3l("see lPoelman & SpaansI (|2005')). It is seen that different lines 
probe different regions in the PDR. The lio — > loi transition 
peaks near the eastern edge of the PDR where the gas temper- 
ature drops towards 50 K, the 3i2 ^ 221 line peaks where the 
temperatures are highest, at the western edge of the cloud. In 
the following we present line profiles when zooming in on dif- 
ferent regions of the PDR, see Fig.|3l The two areas, L and R, 
are indicated by a circle. 

In Figures |5] & we plot the line profiles for a static PDR 
model, and a PDR model with an outflow velocity field, respec- 
tively. Fig. |5] shows the line profiles when S140 is considered 
as a spherical, static cloud. A homogenous (dotted) and inho- 
mogenous (solid) distribution have been considered. The left 
(middle) panel represents the line profiles when one takes the 
mean over a region, with a size comparable to the HIFI beam', 
toward the south-west (north-east) side of the cloud. It is seen 
that different regions lead to different line shapes. In the left 
part of the cloud the abundance of water is lower, due to photo- 
dissociation, compared to the other side of the cloud. This leads 
to a lower water column, hence lower opacity when compared 
to the right side of the cloud. Therefore, a single-peaked line 
profile arises in the left part, whereas the line profile in the 
right part is self-absorbed due to opacity effects. In the inho- 
mogeneous case, the clumps are the opacity source, increasing 
the opacity in the left part, resulting in a top-flatted profile. The 



right panel of Figure|5lgives the line profiles when we take the 
mean over the PDR in its totallity. The shape of the curves is 
similar in both cases. It is seen that the inhomogeneous distri- 
bution produces a higher intensity than the homogeneous dis- 
tribution. 

In Fig.|6lwe plot the line profiles for the 0-H2O 1 10 — » loi and 
3()3 — > 2i2 transitions arising from the static embedded source 
described in Section |2] One notices the difference with Fig. 
|5j in that the larger opacity in the lines causes the profiles to 
be very strongly self-absorbed. Inclusion of a velocity gradient 
with magnitude ~ 3 Avgs does not remove this self-absorption 
characteristic. 

Fig.0shows the line profiles for a PDR model with an outflow 
velocity field. We consider a velocity field, Voutflow ranging 
from -2 to +2 km s"', yielding a velocity equal to zero in the 
middle of the cloud, and gradually increasing towards the edge 
of the cloud, i.e., Voutflow ~ R- As in Fig.|5] a homogenous (dot- 
ted) and inhomogenous (solid) density distribution have been 
considered. The left (middle) panel represents the line profiles 
when one takes the mean over a region, with a size comparable 
to the HIFI beam, toward the south-west edge (north-east) side 
of the cloud. Line profiles coming from the right part have a 
larger width than the ones from the left part. The overall ten- 
dency of the shapes of the line profiles in the inhomogeneous 
case is comparable to the homogeneous models. At scales as 
large as the PDR, (right panel of Fig.0 line profiles are plotted 
in case of a homogeneous (dotted) and inhomogeneous (solid) 
PDR model. Also, a line profile for a model with a smaller 
clump size 1^ = 0.02 pc (dashed-dotted), i.e., with clumps that 
are about 30 % smaller, is plotted. Thi s clump size is consis- 
tent w ith other PDR characteristics (.Spaans & van Dishoeckl 
I1997I) of S 140. We see that the line profile for the homogeneous 
model is more self-absorbed than when an inhomogeneous dis- 
tribution is adopted. Reducing the clump size by a factor of 1.5, 
leads to a more single-peaked line profile, comp arable with ob- 
servations (Snefl et al. 2000:'Ashb v et alb oOOV Although not 
shown, altering the velocity dispersion to 2 km s ' broadens 
the line profile by a factor of ~ 2, in disagreement with obser- 
vations. 

If only microturbulence is considered within the cloud, then 
stronger self-absorption effects occur than when macroturbu- 
lent motions are considered as well. The widths of the profiles 
are broader in the eastern part of the cloud, compared to the 
western side. However, the peak intensity is higher in the west 
side of the cloud. Note that the slopes of the line profiles depend 
on the adopted velocity dispersion. The differences between the 
slopes of the homogeneous and inhomogeneous models do not 
only depend on the density distribution (dumpiness versus ho- 
mogeneity), but are also a result of the chemical and temper- 
ature balance within the clouds. Velocity fields with a moder- 
ate velocity gradient, as in Voutflow, give rise to an almost flat 
line profile with a FWHM of ~ 6 km s~', in agr eement with 
observations dSnell et alJl200 (l lAshhv et al.ll200nh . Similar re- 
sults are found for an infall velocity field with a comparable 
gradient. 



' Half Power Beamwidth between 13" and 39", depending on fre- 
quency band. 
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4. Discussion and Conclusions 

lAshbv et alJ tOOd) presented Monte Carlo radiative transfer 
models for the 557 GHz Ijo — » Iqi ground-state transition of 
ortho-HaO in the S140 PDR observed by SWAS. These mod- 
els adopt power-law representations for their density and tem- 
perature profiles, with the mean density of H2 being a factor 
of ~ 40 larger than the density we adopted, <nH2> = 2 x 10^ 
cm^^, in our PDR models. They are unable to produce a single- 
peaked line profile with the observed integrated intensity by 
any combination of the ortho-HaO abundance and turbulent ve- 
locity when a static cloud is adopted. A combination of turbu- 
lent and macroscopic motions, infall and/or outflow, is needed 
to reproduce the observed line profile. With our adopted den- 
sity and temperature distributions, a single-peaked line profile 
is also not found for a static cloud. When an outflow (or in- 
fall) velocity field is incorporated, line profiles with a FWHM 
of ~ 6 km s ' are obtained for a velocity field with a gradient 
of a few km s However, when modelling the line profiles for 
a static embedded source, with a size ~ 7 times smaller than 
the PDR, i.e., strongly beam-diluted in the SWAS 3.3' x 4.5' 
beam, we find strongly self-absorbed line profiles that remain 
self-absorbed for modest velocity gradients. 
On scales of the PDR, we see that the shapes of the line profiles 
in the inhomogeneous case are less self-absorbed than when an 
homogeneous model is adopted. At smaller scales, and compa- 
rable to the size of the HIFI beam, different line profile shapes 
are found when zooming in on different areas of the PDR. Note 
that the beam of HIFI, depending on the frequency of the tran- 
sition, will be ~ 5 times smaller than the PDR, being able to 
zoom in on certain regions, and to confirm the predicted change 
in the line profiles. Therefore, HIFI can distinguish the embed- 
ded sources from the PDR. The former produce a completely 
different emission line distribution (see Table 3), with the bulk 
of the emission emanating in lines like, e.g., 32i -22i, 3o3 -2i2, 
2 1 2 - 1 1 and 220 - 1 1 1 . Clearly, a combination of pure PDR and 
star-forming core characteristics is needed to plan for and in- 
terpret HIFI data. 



Appendix A: 



for ortho-water abundances, relative to H2, ranging from 10"'" 
to 10"^. The parameters used for this model are listed in Table 
lA.ll Level populations are calculated ona81x81x81 grid 
by means of an escape probability method as described first 
inJPoelman & Spaans (2005). Following Poelman & Soa ansl 
(120051) . we take the distribution of directions, ^, as ^ = 
{-x,x,-y,y,-z,z} for a 6-way partitioned sphere (rather than 



k = {r)) and /3,/x,y,z,ki) = 



for directional labels 



6Tij(x,y,z,ki) 

ki. Comparison is made with the workshop results, and our 
multi-zone escape probability computations are completely 
consistent with the analytical solutions. 
The resulting line profiles for 0-H2O abundances of 10 
10"^, 10 and 10"^ are shown in Figure^] The profiles have 
been calculated for ID (401 grid, slab), 2D (201 x 201 grid, 
cylinder) and 3D (81 x 81x 81 grid, sphere) configurations 
through ray-tracing. The line profiles for the two lowest 
abundances, where the transition is optically thin, are single 
peaked, whereas at high abundances the line profile is no 
longer single peaked but self-absorbed, due to optical depth 
effects. One notices in Fig. |4] that a ID approximation leads 
to absorption effects, in case of a high abundance, that are a 
factor of a few larger than in the case of a 2D approximation, 
which in turn has bigger absorptions in comparison with a 
3D ca se, in agreement with the results described in ljuvelal 
(Il997h . In case of a 2D (ID) approximation, one forces the 
photons to propagate in a disc {along a line), and thus one 
decreases the chance for a photon to escape, leading to larger 
opacities, hence larger absorption effects. The 3D case is in 
good agreement with the analytical benchmark solutions for 
abundances ranging from 10"'° to 10"^ 
In our 'Long Characteristics Method' code, photon packages 

Table A.l. Benchmark model parameters 
Model 



Outer radius 
Inner radius 
H2 density 
Temperature gas 

Gaussian line profile with FWHM 
De-excitation rate coefficient 
Excitation rate coefficient 
No dust or velocity gradient 



0.1 pc 
0.001 pc 
10-* cm-^^ 

40 K 
1 km s"' 
2.18 X 10-'" cm^ 
1.12 X 10-'" cm^ 



The goal of this appendix is twofold. First, we test the escape 
probability code (Poelman & Soaans 2005,) used in this paper 
Second, we present a code based on a 'Long Characteristics 
Method' and compare this code with the escape probability 
approach. 

In order to further te st our escape probability code 
('Poelman & Spaans' '2005"), we try to reproduce the level 
populations and line profiles presented at the water radiative 
transfer workshop, held in Leiden (2004). Three test problems 
were presented: a static two-level problem, an expanding two- 
level problem, and an AGB star with 45 -level ortho- and para- ' ' 
H2O. The first two problems have an analytical solution, see ni(x,y,z)^Rij(x,y,z) = Y^nj(x,y,z)Rji(x,y,z) 
lhttp://www.mpifr-bonn.inpg.de/staff/fvandertak/H20/radxfrtest.pdJ 
and are therefore good test problems, the third is a more real 
istic view. Here, focus is on the first problem. 



travel along fixed directions on a Cartesian grid. 
Consider an atom or molecule having / levels, with Einstein 
Aij, Bij coefficients and collisional rates Cij between levels / 
and j. Let n,(x, y, z) be the population density of the /* level in 
point {x,y,z)- The equations of statistical equilibrium can then 
be written as: 



J*' 



J*' 



where 



We consider the simple case of a two-level water molecule 
(lot and loi states) in an isothermal, constant density sphere 



Rij(x,y,z) = 



I Aij +BijJp + Cijix,y,z) (i>j) 
[BijJp + Cij(x,y,z) (i <j) 



(A.l) 



(A.2) 
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Jp is the total integrated radiation field at position ^(x,y,z) in the 

oo 

cloud, i.e., Jp - J Jp(f>{v)dv. The line profile function, 0,j(v), 



determines at which frequency, i.e., position in the line pro- 
file, the line emits or absorbs, and is normalized to unity, i.e., 

oo 

J (p(y)dv - 1. When thermal as well as turbulent broadening 


are involved, a Gaussian line profile is produced 

1 "'-"o'f 

^ij(y)^- -e^, (A.3) 

with 

Av, = ^(^+^^)'/^ (A.4) 
c m 

where ^ is a root mean-square measure of the turbulent veloc- 
ities. Note that the frequency shift Av = v - vq that enters the 
line profile depends on the velocity field, Av = Av[v(p),v(p')]- 
The angle integrated radiation field at frequency v at position 
^(x,y,z) in the cloud is 

k 

where the sum is over all directions k in the cloud. The number 
of directions is arbitrary, but a 6-ray (4-ray/2-ray) approxima- 
tion is implemented for the 3D {2D/1D) models. The total radi- 
ation field at frequency v at position p(x,y,z) along a direction 
k in the cloud, i.e., ^.1'°' ,, is the sum of the local contribution 
and the contributions coming from the other directions towards 
the point ^(x,y,z)- In case of a ID approximation the upward 
(downward) direction along k is denoted with + (-), i.e., 

p ^+ p" ^- 

Here, t^^, is the optical depth of the line at frequency v between 

two points ^(x,y,z) and ^ (x ,y ,z). The opacity in the line, as- 
sociated with the line profile, can be written as 

a,V = ^rijix, y, z)Bji(l - ^)cf>ij(v). (A.7) 
An ' njgi 

The optical depth follows by multiplying the opacity with the 
distance dl = p(x,y,z) - p (x ,y ,z ) between the two gridpoints. 
The local intensity initially carried by a photon package in a 
point p(x,y,z) is given by 

''^'p%.y,z) = ^Aijhvijmix, y, z)4>ij{v)dl + Dust + CMB, (A.8) 

with dl the size of the corresponding gridcell. 
In this way, we use global information to determine iteratively 
in every gridpoint the ambient level populations. When test- 
ing this code against the benchmark model described in the 
beginning of this appendix we find that the level populations 
are completely consistent with the analytical solutions, as do 
the level populations calculated with our escape probability 
method. However, the long characteristics scheme, even though 
it always works, is very time consuming when a large grid is 
involved. 



Note that the escape probability method used in this paper, as 
described in Poelman & Spaans (2005) is not a one-zone but a 
global approach. It has been modified such that the calculation 
of the level populations in a grid point depends on level popula- 
tions in all other grid points through the different escape prob- 
abilities that connect any two grid points. This multi-zone im- 
plementation of the escape probability constitutes a large accel- 
eration over long characteristics and Monte Carlo techniques. 
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Velocity (cm s ^) Velocity (cm s ') Velocity (cm s ^) 

Fig. 4. Line profiles for a static two-level problem (Iqi and loi states). The line profiles are calculated for a ID {left), 2D {middle) 
and 3D {right) configuration for different 0-H2O abundances. The bottom curve resembles the line profile when the 0-H2O abun- 
dance is 10"^*^, whereas the top one is the Une profile in case the 0-H2O abundance is 10"'. Y-axis in units of erg cm"^ s"' sr"' 
per cm s"' . 




Fig. 5. Line profiles for the 1 10 —> loi transition in case of a static PDR are plotted. Line profiles are plotted when zooming in on 
a certain area in the PDR (left and middle panel), and when taking the mean over the PDR (right panel). The left {middle) panel 
represents the south-west {north-east) part of the cloud. Dotted {solid) lines are the line profiles in case of an (m)homogeneous 
model. Y-axis in units of erg cm"^ s"' sr"' per cm s"^ 




Fig. 7. Line profiles for the 1 10 ^ loi transition are plotted when a systematic outflow velocity field is employed. Line profiles 
are plotted when zooming in on a certain area in the PDR (left and middle panel), and when taking the mean over the PDR 
(right panel). The left {middle) panel represents the south-west {north-east) part of the PDR. Dotted {solid) lines are the line 
profiles in case of an (i>i)homogeneous model. In the right panel a fine profile for a model with a smaller clump size, Ic = 0.02 
pc, (dashed-dotted) is also plotted. Y-axis in units of erg cm"^ s"^ sr~^ per cm s"^ 



